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Abstract A variety of phenomena in physical and biological sciences can be math- 
ematically understood by considering the statistical properties of level crossings of 
random Gaussian processes. Notably, a growing number of these phenomena demand 
a consideration of correlated level crossings emerging from multiple correlated pro- 
cesses. While many theoretical results have been obtained in the last decades for in- 
dividual Gaussian level-crossing processes, few results are available for multivariate, 
jointly correlated threshold crossings. Here, we address bivariate upward crossing 
processes and derive the corresponding bivariate Central Limit Theorem as well as 
provide closed-form expressions for their joint level-crossing correlations. 

1 Introduction 

Various phenomena in the biological or physical sciences are amenable to the de- 
scription by level crossings of random Gaussian processes [1,2]. Examples of these 
phenomena are spike coordination of neurons in the brain [3], insurance risk assess- 
ment [4] and stress levels generated by ocean waves [5]. Therefore a number of math- 
ematical studies in recent decades have focused on the statistical properties of level 
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crossings arising from stationary Gaussian processes [2]. However, largely this liter- 
ature addresses the properties of one level-crossing process and rarely deals with the 
coordinated level crossings of multivariate Gaussian processes. A prominent appli- 
cation where correlated level crossings are of particular importance is neuroscience. 
Recent work has shown that the spikes of a cortical neuron can be approximated by 
a Gaussian level-crossing process [3, 6]. The assumption of Gaussianity is prompted 
by the experimental observation that cortical neurons are on average connected to 
~ 10000 neurons and therefore receive a barrage of inputs that together lead to a near- 
Gaussian fluctuation at the cell body of any given cortical neuron [7]. The spikes of 
two neurons are then modeled as upward level-crossing times of two cross-correlated 
fluctuating Gaussian potentials. 

In this article we aim to address two features of level crossings of multiple cor- 
related Gaussian processes. First, we want to clarify whether level-crossing counts 
derived from multiple correlated processes are jointly Gaussian. Second, we want to 
understand how many more coincident level crossings in a given time instance are 
expected if the underlying Gaussian random processes are correlated. Let us provide 
an intuitive reason for these questions. Starting with the first question, we recognize 
that if level-crossing counts of two neurons were jointly Gaussian, then a simple mea- 
sure of dependence is the covariance or the Pearson correlation coefficient. Measur- 
ing a vanishing correlation coefficient or vanishing covariance between two neuronal 
spike counts would in this case imply true statistical independence, because only in 
the case of multivariate Gaussian distribution is it permissible to conclude indepen- 
dence from vanishing count correlations. This implication is not permissible if the 
marginal distributions are not Gaussian or are Gaussian but the joint distribution is 
not a multivariate Gaussian distribution. While marginal Gaussianity has been shown 
for level-crossing counts in [2] for large bin sizes, joint Gaussianity is still an open 
question. It might seem natural to imply joint Gaussianity from marginal Gaussianity 
for multivariate level-crossing processes, however, numerous counter examples exist 
to prove this intuition wrong, see Sect. 5 in [8]. Here, we use a modified Breuer- 
Major Theorem to prove joint Gaussianity and show that any linear combination of 
level-crossing counts of the two processes is also Gaussian. 

The second question we address in this article deals with the conditional proba- 
bilities of two level-crossing processes. We are interested in how the level crossings 
of one Gaussian process can be used to predict the level-crossing probability of the 
partner process in a specific time interval relative to the observed level crossing in 
one process. In neuroscience, coordinated neuronal firing drives changes in synaptic 
connectivity and calculating the spike count dependencies across neurons is therefore 
a topic of current research efforts (e.g. Chap. 8 in [9]). The available mathematical 
results for conditional upward crossings in Gaussian processes currently comprise 
mostly variance and moments for one level-crossing process (see Chaps. 3-5 in [2]) 
as well as the low and high correlation limit in pairs of processes [3, 10]. As yet, 
a comprehensive closed-form solution covering the complete level-crossing cross- 
correlation function is currently lacking. Here, we use a regression approach to de- 
rive, for all correlation strengths, the conditional level-crossing correlation functions 
in two continuous Gaussian processes. We hypothesize that the level-crossing cor- 
relations we provide in this article could also be valuable in other fields outside of 
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Fig. 1 Cross-correlations in the Gaussian variables lead to correlations of coincident level crossings. 

a Spike correlations can arise from common input in a neuronal network, b We consider coincident level 

crossings arising from two Gaussian processes that share a common latent source. Whenever the voltage 

crosses a threshold \J/ from below a spike is emitted. Spikes are indicated by vertical solid lines. The 

y. 

vertical dotted lines indicate the width of a time bin T used to compute spike counts Uvq ^ (ifr), i = 1,2 



neuroscience for example in risk assessment calculations to predict the risk of joint 
default for insurance purposes. 

The article is structured as follows. In Sect. 2 we define the mathematical model 
setting and introduce the concept of level crossings and specifically the upward cross- 
ings. In Sect. 3 we use a regression approach to obtain a general closed-form solu- 
tion for cross-correlations of level crossings in two correlated Gaussian processes. 
In Sect. 4 we prove the joint Gaussianity (Central Limit Theorem) for the correlated 
joint upward crossings for two correlated Gaussian processes. In the section on mate- 
rials and methods (Sect. 6) we provide detailed derivations of the reported results. We 
assume throughout this article that both level-crossing processes arise from crossings 
of the same threshold level by two Gaussian processes with different variances. This 
is permissible because the number of level crossings, the Rice rate [11], depends only 
on the variance-to-threshold ratio, but not on these quantities individually. We there- 
fore work with a pair of level-crossing processes where each process has a unique 
voltage variance and therefore the rate of crossings in the two neurons being con- 
sidered are, unless stated otherwise, not the same. Let us note that this assumption 
is prompted by the observation that in a living brain typically no two neurons are 
identical in all their properties and differ at least in their firing rate. 



2 Mathematical Definitions of Multivariate Level Crossings 

Here we address the statistics of coincident level crossings arising from two Gaussian 
processes that share a common latent source. This situation is illustrated in Fig. 1(a). 
We choose to illustrate the situation using a neuroscience perspective. Neurons in a 
mammalian brain receive synaptic inputs, both excitatory and inhibitory, from thou- 
sands of other neurons. Particularly in the visual cortex, the excitatory and inhibitory 
inputs largely cancel each other and lead to a net fluctuating residual current at the 
cell body of each neuron. These residual fluctuations drive the spikes of individual 
neurons. These voltage fluctuations arise from largely independent inputs so they are 
well approximated by a random Gaussian process with temporal correlations deter- 
mined by the temporal structure of synaptic currents [7] . 
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2. 1 Definitions of Multivariate Voltage Distributions 

We begin by defining the random, temporally correlated Gaussian zero mean process 
Vj(t) which represents the voltage of a neuron j 

/oo 
e ia fl' 2 (A)(Vr^7d^(X) + VFdW c (A)), (1) 
-00 

where fy is a combination of filters fy (A) = y (A)g(A), where y represents the mem- 
brane filter and g the synaptic filter. Both of these filters can be chosen freely, but 
their combination should guarantee a continuously differentiable voltage trajectory. 
Wj with j = { 1 , 2} and W c are complex random measures with independent incre- 
ments, such that for all Borelian sets A e S(R) we have E|W(A)| 2 = m(A), the 
Lebesgue measure of A. By W c we denote the common noise component. Moreover 
if A n B = 0 then W(A) and are independent Gaussian random variables. The 

correlation strength r, r = [0, 1), denotes the presynaptic overlap of neurons 1 and 2 
generated in a neuronal network, and it is illustrated in Fig. 1. If r = 0 the voltages 
V\ and V2 are independent if r = 1 the voltages V\ and V2 are identical. The auto- 
and cross-correlation functions between V; and Vj are, respectively 

CvjW = {Vj(0)Vj(T)) = or* c(r), (2) 

= (V/(0) V*(r)) = rory.or y ,c(r) for j, k e {1, 2}, (3) 

where c(r) = e lxX fy (A)dA,, and r is the considered delay. The vector (Vi(0), 
7/(0), V2(r), V^W) comprising the voltages and their derivatives is Gaussian and 
has the covariance matrix 

0 Z 13 (t) S u (r)\ 

1 9 , (4) 

where the variances are cry. = Cy- (0), oy, = — Cy. (0) and covariance functions are 

j j j 

given by 

Sl3(T) = (Vi(0) V 2 (T)) = rcr yi cry 2 c(T), 

A 4 (r) = (Vi (0) V 2 '(t)) = rcr^cry.c^r), 
^24(r) = (V/C0) V 2 7 (t)) = -ro Vl o Vl c"(x). 
We use the correlation time r s to quantify the width of the correlation function c(r): 

x s = ^c{x)/\c"{x)\. (5) 

If the filters y(X) and are classic low-pass filters, then the correlated voltage 
processes of Eq. (1) can be written in a differential form for each neuron j: 

TMV}(t) = -Vj(t) + Ij(t), (6) 



( < 

0 

^13 (r) 
^I4(r) 
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where Ij (t) is the residual Gaussian current fluctuation with variance a, xm the mem- 
brane time constant of the neuron, e.g. [12-15]. The synaptic drive Ij(t) can be sep- 
arated into two parts: a common and an individual noise component 

/oo 
e ifA yo(X)(yr^7dWi w + ^&w c {\)) , a) 
-oo 

where Q (A) is the synaptic noise spectral density. Using Eq. (6) we obtain the fol- 
lowing spectral representation for the stationary solutions: 

/OO lQ Q\ 

.oo (l + /r M A) v ' 

In this form the spectral density of each V) is given by fy(X) = Q(A)/(1 + r^A 2 ). 
Analogously to Eq. (3), we obtain 

/OO Q fy\ 

e irk - dk = ra l a 2 c(z). (8) 

-oo (1 + T M A ) 

2.2 Upward Crossing Definitions 

Neurons communicate using brief pulses, the so called spikes, which are emitted 
whenever a voltage threshold is crossed [9]. The integrate- and- fire- type neuron mod- 
els that are frequently used in computational neuroscience [12-15] generate a spike 
in neuron j any time a voltage Vj(t) crosses a fixed threshold ifr and subsequently 
reset the voltage to a reset potential. Recently, it has been shown that in many phys- 
iologically relevant cases the leaky integrate-and-fire model can be equivalent to a 
level-crossing model without reset, where spikes are modeled as positive threshold 
crossings and are not followed by a reset [3, 6, 16]. Here, we therefore identify the 
spikes of a neuron j with the positive level crossings of its voltage Vj(t) and quan- 
tify the cross-correlation between level crossings in neurons 1 and 2 by the following 
level functional: 

Uq( B )W) = ] im n 7^2 f U\ViW-f\<8}Vl(si)l {V r {s) > 0] 

8^0 (Ad)^ jQ( e ) 1 

x 1 {|V2(^l+^2)-^l<5} V 2^1 + s 2)l{V^ Sl +s 2 )>0}dsidS2, (9) 

where Q(e) = I\ x I2 is a bounded and finite rectangle in R 2 , \jr denotes the voltage 
threshold in both neurons (see Fig. 1). Here 8 is introduced to quantify the infinites- 
imal interval around the threshold \J/ where a spike takes place. We choose the same 
threshold for both neurons and two different variances (ay l ^ ay 2 ) and keep all other 
parameters the same. oy x ^ ay 2 represents the biological situation in which two neu- 
rons of the same neuronal type could have differences in the strength of their synaptic 
input and threshold-to-variance ratio but are exposed to the same temporal back- 
ground statistics. We will consider the following random field Z : R 2 x Q — ► R 2 , 
defined as (s\,S2) -> Z(s\,S2) = (Vi(s\), V2(s\ + S2)). & denotes the probability 
space; here Q is the Gaussian probability space. The field Z is Gaussian and Z(s\ , si) 
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and Z(0, S2 — ^1) are equal in distribution. We denoted by p S2 -s\ (•> ■) the bivariate 
Gaussian density of vector (V\(s\), V^fe))- If Q(s) = [t 9 t + s]x [z, r + s] and the 
prerequisites of Theorem 6.2 in [2] are fulfilled we can write 

E[t/ G( *)«0] 

= e G(e) : V1O1) = ^, V 2 0i + * 2 ) = ^ 



E[ I detZ 7 (si,s 2 )\, 1v'(ji)>o. Iv^i+^olZfai.^) = VO] 
x p S2 - Sl (f, ^)d^id^ 2 (11) 

= s E[|det Z'(0, s 2 ) I ly/ ( o)>o. 1 v 2 '(, 2 )>ol z (0> J 2 ) = VO] 

xp, 2 (^^)d5 2 (12) 

= ^ r+£ E[|y 1 / (0)||y 2 / fe)|ly; (0 )>o, lv 2 '(* 2 )>ol^i(0) = Viisi) = if] 

x p S2 (if, i/r)ds2, (13) 

where the expectation value is denoted by E, and det(Z(,si, s^)) is the determinant 
of the correlation matrix for the vector field Z (s\, S2). Now, we are left to prove 
the conditions of Theorem 6.2 in [2]. First, we find that conditions (i) and (ii) of 
Theorem 6.2 are satisfied by definition. Condition (iii) holds because p S2 (if, if) is 
not degenerate. If we let I\ and 7 2 be two finite and bounded intervals in R, condition 
(iv) is satisfied because 

P{3(ji, j 2 ) 6/ix/ 2 : Z(s u s 2 ) = (if, if), detZCsi, s 2 ) = 0} 

< P{*i g /1 : Vi(si) = ^r, 7/(51) = 0} + P{* 2 e h : V 2 (s 2 ) = if, V^(s 2 ) = 0}. 

Here P denotes the probability measure. We can define the correlation of two spike 
trains as 

(si(t)s 2 (t + r) := lim 

X Pr (if,if) (14) 
VlV2Pr(if, V\, if, \)2)&V\&V2, (15) 



/•OO /*C 

Jo Jo 



where /^(^ ^ ^2) is the joint Gaussian density of the vector (Vi(0), 7/(0), 
V2CO, V^r)). The conditional firing rate v coik i(t) then is 

Vcond(^) = [s\(t)s 2 (t + r))/Vviv 2 , (16) 
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where vj = 2jr J v exp(— ^r) is the firing rate of a neuron j, for j = 1, 2. In the next 
sections we provide closed-form expressions for {s\(t)s2(t + r)> and v con d(r). 



3 Cross-Correlations of Two Upward Level Crossings 

Here, we address (s\ (t)s2(t + r)) and provide a closed-form solution that is valid for 
any cross-correlation strength r between two level-crossing processes and any time 
delay r. 

Proposition 3.1 Following the steps outlined in the methods section, Sect. 6.1, we 
can apply a regression model and Mehler's Formula (see Lemma 10.7 in [2]) to prove 
that 

(si(t)s 2 (t + r)) = C(a ib) (r)p r (1r, if), (17) 

where p T (-, •) is a joint Gaussian distribution for voltages Vj as defined in Sect. 2.2 
and C(a,b)( x ) is the series given by 

C (a ,b)M = [b#(b)-4>(b)][a0(a)-4>^ 

^ n\(a 6l (r)a €2 (r)r- 1 

where 

<r 6l (T) = (<r£, - +2aia 2 i:i3(r) +a^ 2 )) 1/2 , (19) 

or 62 (r) = (or| - (/^ + IfahZnW + ^ 2 2 ^ 2 )) 1/2 , (20) 
Cov (€lf€2 )(T) = 172400 - (ftaior^ + (aiA2 + a2^l)^i3(r) +a2fta^y 2 ), (21) 

-(^(ai+a 2 )) , -(^(fr+jfe)) 



cr 6l (r) a 62 (r) 



27i 4 (T)i:i3(T) -Uu(r)cr^ 



d\ = — ^ — 5 — — - r, «2 = 



or^or^ - i7i3(r) 2 ' a 2 ^ - ^(r) 2 ' 

<4 2 ^i4(r) n -U 14 (r)E 13 (r) 



Pi = — 9 ~ — TT^"' ^2 



a 2 ^ - ^(r) 2 a 2 ^ - I7i 3 (r) 2 

0 0 are the standard Gaussian density and distribution, respectively, and 
0 = 1 - 0. H n (z) = (~^) n -^( e ~ z /2 ) eZ 11 are tne Hermite polynomials. Note that 
the first two terms in Eq. (18) correspond to truncation orders n = 0 and n = 1, 
respectively. 
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Fig. 2 Convergence ofpairwise level crossing correlations, a Spike correlation function v CO nd( r ) vs - time 
lag r for different truncation orders n in Eq. (17); r = 0.7, cry l = cry 2 = 10 mV, r s = 20 ms, v = 5 Hz, 
^ = 9.64 mV. b v conc i(T) vs. x for a pair of rate heterogeneous neurons with r = 0.7, oy l = 10 mV, 
ay 2 = 5 mV, r s = 20 ms, v\ = 5 Hz, V2 = 1.24 Hz, ^ = 9.64 mV. c v conc i(T) vs. time lag for varying 
correlation strengths r in a pair of neurons with oy x = ay 2 = 10 mV, x s = 20 ms, v = 5 Hz, ^ = 9.64. 
d v con( }(T) vs. time lag for varying correlation strengths r, in a pair of rate heterogeneous neurons with 
ay l = 10 mV, oy 2 = 5 mV, r s = 20 ms, v\ = 5 Hz, = 1.24 Hz, ^ = 9.64. The truncation order of 
v con( i(r) -series in Eq. (17) in c-d is n = 10. In a-d the filled circles at t = 0 indicate the predicted 
v cond(0) ( as m Eq. (23)) and colored squares denote the corresponding numerical simulations obtained 
with TV = 2000 independent realizations of 20 s length 



To aid the explicit evaluation of C( fl ^)(r) in Eq. (18) we provide code for a com- 
puter algebra system. 1 Figure 2(a), (b) demonstrate v con d(r) obtained using Eq. (17) 
for progressively large truncation orders n. Notably, we find close correspondence 
between the first truncation order n = 1 and the large n limit (n = 10). Figure 2(c), 
(d) show Vcond(^) vs. r as in Eq. (17) for varying correlation strength r. For sim- 
plicity, we chose c(r) = cosh(r/r 5 ) _1 and r e [0, 1). We note that v con d(r) for two 
identical neurons (cr^ = ay 2 ) is a symmetric function while for a pair of neurons with 
different rates (ay 1 ^ cry 2 ), v C ond(^) is asymmetric. 

Let us now briefly discuss the result we obtained in Eq. (17) within the context 
of previous level-crossing literature. One of the closely related results is the vari- 
ance of level crossings and maxima provided in Proposition 4.4 in [2]. However, this 



We provide the MATHEMATICA 8 (Wolfram Research) code to iteratively calculate v con( }(T). The code 
can be found at: www.tchumatchenko.de/CodeNuCond_Fig2.nb. 
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result is derived for one level-crossing process, while we addressed a pair of level- 
crossing processes. For multiple cross -correlated processes orthant probabilities that 
describe expressions of the form F(V\ (t) > ty\ , V2W > ^2) have been obtained, e.g. 
Lemma 4.3 on p. 78 in [2]. However, specific results for the cross-correlations of 
upward crossings are sparse. For two correlated upward level-crossing processes pre- 
vious studies have addressed the limiting cases of weak (r ~ 0) or strong (r ~ 1, 
r < 1) cross-correlations [3, 10]. However, to address upcrossing correlations in the 
intermediate regimes where neither r ~ 0 nor r ~ 1 no analytical methods are avail- 
able. Therefore, it was previously necessary to numerically evaluate the associated 
Gaussian probability densities in Eq. (15). The direct numerical evaluation of multi- 
dimensional Gaussian integrals can be computationally and algorithmically demand- 
ing, requires adaptive grid procedures and its accuracy can be hard to evaluate [17]. 
For the specific case of r = 0 we show in the materials and methods section on 
'Zero time lag correlations', Sect. 6.2, that a direct evaluation of the Gaussian double 
integral is possible via a variable substitution. The key to this variable substitution 
method was a manageable unity correlation matrix. For any other finite r > 0 and a 
given finite correlation strength 0 r < 1 we could not identify a transformation that 
leads to a tractable integral and we therefore derived the series expansion in Eq. (17). 
This solution is explicit such that each series term of order n can be evaluated and 
studied separately. Furthermore, Eq. (17) consists of analytical functions with a well- 
studied parameter dependence. This makes it possible to predict the influence of a 
specific parameter, such as time scale r s , firing rate v\ or correlation strength r on the 
upward level-cross correlations. As an example, we evaluate Eq. (18) for an identical 
pair of neurons up to the third order in r via a Taylor expansion. We obtain 



We recognize that the linear and quadratic expressions are equivalent to the first 
and second order r -expansion reported in [3, 10]. The cubic term has not been re- 
ported before, to the best of our knowledge. This demonstrates consistency with pre- 
vious results and illustrates that expansions of any order can be obtained via Eq. (18). 
In this context, it is desirable to have an exact reference point for deciding how many 
n- orders are necessary for Eq. (17) to be sufficiently accurate. Such a reference point 
can be the zero lag value which we calculate exactly. Deviation of Eq. (17) for a spe- 
cific n from this reference point can serve as an accuracy benchmark. Following the 
calculations in the methods Sect. 6.2, we derive 



VcondM = v + rv(c(x)ir 2 /al - itx^c' ' (x) /2) 





2 




(22) 
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Vcond(O) = 



1 + (2r arctan(V(l + - r)))/Vl r 2 



47r 2 v /viv 2 r 2 




(23) 



Figure 2(a), (b) demonstrates v C ond(^) obtained using Eq. (17) for different trunca- 
tion orders n alongside the zero lag correlation v con d(0). Figure 2(c), (d) demonstrates 
v con d(T) obtained using Eq. (17) as a function of the correlation strength r alongside 
the zero lag correlation v conc i(0). As previously, we chose c(x) = cosh(r/r^) _1 and 
r G [0, 1). We note that for two identical neurons (cry l = oy 2 ) v con d(^) is a symmetric 
function. Yet, for a pair of neurons with different rates (oy x ^ oy 2 ) the spike correla- 
tion function v con d(r) is asymmetric, indicating that the lower rate neuron spikes on 
average after the higher rate neuron. 

3.1 Relation to the Leaky Integrate- and-Fire Model 

Here, we address the relation between spike statistics and spike correlations in the 
level-crossing setting in our article and previous results in the leaky integrate-and- 
fire framework [12-15]. In a current-based leaky integrate-and-fire model driven by 
fluctuating noise the voltage of a neuron Vj (t) is given by 



where Ij(t) is the input current of a neuron, ^(t) sl white noise, unit variance drive. 
The voltage power spectrum for this model is a combination of low-pass filters 
fy(X) ~ [(1 + r^A 2 )(l + r 2 A 2 )] -1 and its correlation function can be determined 
according to Eqs. (8). If the voltage Vj reaches the threshold 0 the neuron j emits 
a spike and the voltage is subsequently reset to a reset value V r . The integrate-and- 
fire model differs only in one important detail from the level-crossing approach — the 
presence of a reset after a spike. A recent article by Laurent Badel systematically 
compared the validity of upward level-crossing approximation for the firing rate, 
spike correlations and frequency response of a leaky integrate-and-fire neuron [16]. 
This study reached the conclusion that the upward level-crossing approach accurately 
represents the leaky integrate-and-fire model if two conditions are fulfilled: (1) the fir- 
ing rate is much lower than the typical relaxation time of the voltage, (2) the synap- 
tic filtering time constants remain of the same order of magnitude as the membrane 
time constant (tj/tm ~ 1). Numerically, the validity of the approximation remained 
highly accurate even for synaptic time constants 0.4<t//tm<2.6. 

A number of spike correlation results have been derived in the leaky integrate-and- 
fire model for the limit of weak correlations [18-20]. They include the observation 
that the spike correlation coefficient increases with firing rate [18, 19]. The equivalent 
firing rate dependent increase in spike correlations and correlation coefficients for 
low correlation strengths has been reported for the level-crossing model, see [3] and 
Fig. 3(A) (right) and Fig. 2(A) (top) in [10]. Furthermore, leaky integrate-and-fire 
model exhibits a sublinear dependence of correlation coefficients on input strength r 
[18, 21], which we see confirmed in Fig. 4. 



T M Vj(t) = -Vj(t) + Ij(t), 

T I I , j (t) = -Ij(t) + a%(t), 



(24) 



(25) 
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(a) (b) (c) 

Marginal Gaussianity joint Gaussianity test (linear combinations) joint Gaussianity test (Mahalanobis) 




U^Y) - E[U£ T] (i|J)] Angle 9 D 2 x 2 



Fig. 3 Univariate and multivariate count Gaussianity. a Univariate Gaussianity of counts T ^{^) — 
EfC/^Q 1 j,j for large bin size T = 25t s . A solid black line represents the corresponding zero mean Gaus- 
sian fit. b Shapiro test p-value of the projections A = xj • (cos(#), sin(#)), for all 9 e [0, 2tt). c Prob- 
ability density of D 2 (left). QQ-plot between the theoretically predicted x|-quantiles and the empirical 
quantiles of D 2 in Eq. (35). Figures b and c are both validations of the multivariate Gaussianity of the 
bivariate vector X in Eq. (34). In all panels \(r = 0.3, x s = 1 ms, ay. = 1 



4 Joint Gaussianity of Upcrossing Counts 

Spike count cross correlations and correlation coefficient measurements in pairs of 
neurons are ubiquitous in neuroscience and are often used to measure the strength of 
interdependencies in a pair of neurons, e.g. in cortical neurons [18, 19, 22], in model 
neurons [23] and in theoretical and experimental studies of net correlations emerging 
in recurrent networks [23-28]. Spike counts and their cross correlations in neuro- 
science are often computed for a variety of bin sizes varying from T = 0.1-1 ms [22] 
to T = 2 s [29] . Here, we are interested in the question when spike count correlations 
of two neurons computed in a bin size T are jointly Gaussian such that their cross 
correlations are unbiased measures of statistical dependence or independence. 

In this section we will prove that the spike counts of two neurons, approximated 
by up crossings of a Gaussian process, approach a joint multivariate Gaussian dis- 
tribution for large bin sizes T. We start by considering the marginal distributions 
of upcrossing counts. From the one-dimensional Central Limit Theorem proven in 

[2] we know that -^=(U^ — E[t/^ r j(V r )]) converges for T —> oo to a one- 
dimensional centered normal variable with finite variance. We provide a direct illus- 
tration of this classical univariate result in Fig. 3(a). Now, it is tempting to conclude 
that because the distribution of counts in each neuron is Gaussian, the joint distribu- 
tion for the vector (U^ T ^ (x//), f/^V] (^)) ^ s a ^ so a multivariate Gaussian distribution. 
Yet, this conclusion is mathematically forbidden. While a joint Gaussian distribution 
implies marginal Gaussianity it is general not possible to inverse this relation [8]. 
The joint Gaussianity of spike counts is a highly desired property. If two counts are 
jointly Gaussian zero count correlation directly implies statistical independence. If 
count correlations between neuron 1 and neuron 2 are zero such that U^ T ^(\jr) and 

U^ T j(\f/) are uncorrected, then only if the vector (U^ T ^), U^ T ^)) is from a 
multivariate Gaussian distribution is it possible to infer independence of neuron 1 
and neuron 2. Let us consider a teaching counter example where vanishing correla- 
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tion between the variables X and Y does not imply independence: X e A/"(0, 1) and 
Y = X 2 . We obtain Cov(X, Y) = 0, but the two random variables are strongly linked. 
Contrasting examples of where X and Y variables are both marginally but not jointly 
Gaussian, have a zero correlation but are not independent can be found in Sect. 5 
in [8]. To benefit from joint Gaussianity and be able to infer true statistical indepen- 
dence from vanishing correlations, we prove here the joint Gaussianity of upward 
level crossings/spike counts for large T . In the following we derive the Central Limit 
Theorem for the spike counts of two neurons, using two steps. First, we use the one- 
dimensional Central Limit Theorem proven in [2]. Second, we apply a version of the 
Breuer-Major Theorem adapted to our upward crossing setting which we present in 
Sect. 6.4. 



Theorem 4.1 Let Vi(t), i e {1,2} be two processes satisfying Eq. (1), with covari- 
ance Cy tj (t) = E[Vj (0) V)(r)] where i, j = 1,2 and a common spiking threshold 
To take advantage of the available Gaussianity proofs that are typically derived for 
unit variance processes, we rescale the voltages Vi(t) and the threshold iff to ob- 
tain processes Xi (t) with unit variance and unit variance of the derivatives. Xi = 



Vi(yc v ..(0)/|Cf..(0)|0 



then has the correlation function c(x) 



Cv I - f (yc v ..(0)/|C , ;..(0)|T) 



and the spiking thresholds are = ^r/ y /Cy ii (0). The number of tyi -lev 'el upcross- 

V. 

ings in a time interval T for process Xi is given by U^ l T ^(^i). We assume that the 

following necessary and sufficient conditions hold. First, K{(U^ l T ^(\l/)) } < oo. This 
holds if and only if c(z) satisfies Geman's Condition (see, e.g., Theorem 3 in [30]). 
Second, / 0 °° |c^(r)|dr < oo where n £ {0, 1, 2} is the order of derivation. Then 



1_ ( U [0,T] 



[0 f V](^l)] 



Vf \u* 2 




(26) 



, [05r] (^2)-E[^ r] (^ 2 )]y 

where the count covariances a\j with i, j £ {1,2} are three convergent series. Each is 
then given by 0 < an = YlqLi °\ i (#) an d 0 < a \2 = Y^q=\ <7 x u x 2 (q)> b° tn of which 
are finite. The first two terms in these series for \jr = \jrt, Cy u (0) = Cy 22 (0) are 



J x 1 



(l) = 24>m 



2 \r 



a 2 (2) = 20W 



+ 



2 

|_2tt 

\_2ix 
1 1 

7T ~ 4 



pOO 

JO 

{f 2 



c(s)ds 



4 Jo 



c"(s)ds 



(27) 



"7 C 

Jo 



f 1 



c(s) 2 ds 



c'(s) 2 ds 



+ 



^ poo ^ poo 

— / c"(s) 2 ds --ir 2 I c(s)c"(s)ds 
2n Jo 4 Jo 



(28) 



(29) 



and the first and second order covariances are 

,2 



^ J c(s)ds --J c"(s)ds , (30) 
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a M (2) = 20 2 (iA)r 2 



i 9 r°° 

-(**-!) Jf c(,) 2 d, 



"CO 

+ {{i-i)+ 2 -i)] 0 c ' isfds (31) 

] /*oo ] poo 

+ — / c"(,y) 2 dy- -f 2 / c(j)c // (j)ds . (32) 

We note that only for r = 1 we obtain <Jxi,x 2 U) = °\ x O)- asymptotic Pearson 
correlation coefficient pj , defined by 

^12 

hm p r = - _J — = (33) 

yVar(^ r] (^i))yVar(^ 2 r] (^2)) 

w/// a/s6> converge to the respective ratio of the asymptotic covariances and variances. 

4.1 Numerical Confirmation of Joint Gaussianity and Limit Covariances a\j 

In the last section we showed that spike counts of two neurons in large bins approach 
a bivariate Gaussian distribution with finite variances. Here, we illustrate this theoret- 
ical result in simulations of level-crossing processes. We choose two methods based 
on linear combinations and the Mahalanobis distance to numerically confirm joint 
Gaussianity. First, we numerically confirm the joint Gaussianity by showing that all 
linear combinations of two simulated spike counts for a large bin are Gaussian. We 
consider a vector 



X:= 



(^«V]W -4<Vi«o]), ^« 2 r]W -EfeiW])) ( 34 > 



where [/^^(i/r), U^ T ^) and ElU^^x//)], E[U^ 2 t ^(\I/)] are the spike counts 
and their average in neuron 1 and 2, respectively. X consists of N x two-dimensional 
samples. N denotes the number of sample realizations and i e [1, N] the consecu- 
tive sample number. We project each two-dimensional element X/ in all directions 
(cos(0), sin(0)) by calculating the scalar product A = Xj • (cos(0), sin(0)), with 
6 e [0, 2tv). Subsequently, we apply a Shapiro test on A to verify Gaussianity for 
all univariate projections. The p -value of this Shapiro test conveys the certainty with 
which the Gaussian hypothesis cannot be rejected. As a second numerical test of 
joint Gaussianity we use the Mahalanobis distance. This test is based on the fact that 
if X ~ A/j(/x, U), where d is the dimensionality, /x is the mean and U the standard 
deviation, then the Mahalanobis distance D 2 with entries Df 

Df = {(Xi - iL)'E-\Xi - /x), * = 1, . . . ,n], (35) 

is distributed according to a x j -distribution with d degrees of freedom (see, e.g., 
Sect. 3.1.4 and Eq. (3.16) in [4]). By numerically estimating the count sample average 
/jl and U we calculate in our case Df and compare it with a xf -distribution, using 
the QQ-plot method (see Fig. 3(c)). 
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10 15 
time bin T 



Fig. 4 Finite asymptotic limit of covariances a^j and correlation coefficient pj . a Count covariance 
a\2 = Yl^Li a X\ X 2 ^ vs ' ^ me ^ m ^ ^ or var yi n § correlation strengths r in the case of two identical 
neurons. The covariances converge for increasing time bin T towards the value predicted in Theorem 4.1 
(indicated by small, thick lines). Here x s = 1, and iff = 0.3. b Limit value of the covariance a \2 vs. r for 
T = 20t s (from a) vs. r. The case of r = 1 (black) corresponds to the variance a\\ ~ 5.5 • 10 -3 , indicated 
by the dashed horizontal line, c Correlation coefficient pj vs. r for large time bins (T = 20r s ). The dashed 
line indicates the equality line 



Figure 3 demonstrates the results of the joint Gaussianity tests for a bin size 
T = 25 r s , where z s = 1 ms. Figure 3(a) shows the empirical univariate distribu- 
tion of spike counts in one level-crossing process derived from N = 10000 inde- 
pendent count realizations. Figure 3(b) demonstrates that in N = 10000 independent 
count realizations of X p -values for all 6 are above the 10 % significance level. Fig- 
ure 3(c) (left) illustrates that the Mahalanobis distance D 2 of a two-dimensional spike 
count variable X are well approximated by the x| -distribution (solid line). Figure 3(c) 
(right) demonstrates in a QQ-plot of the empirically measured Z) 2 -quantiles and the 
theoretical X2 -quantiles that they are linearly related. This is an indication that both 
distributions are equal. 

Figure 4 addresses the numerical confirmation of the constant asymptotic covari- 
ances aij introduced in Theorem 4.1 in Eq. (26). We choose \jr = 0.3, x s = 1, for 
N = 5000 independent count realizations. To numerically compute the covariances 
aij from spike count simulations we used the covariance-matrix estimator proposed 
by [31]. Figure 4(a) demonstrates the convergence of covariance a\2 to a finite value 
that is predicted by Theorem 4.1. The asymptotic large T limit is denoted by a brief 
colored horizontal line. Figure 4(b) demonstrates the dependence of this asymptotic 
limit on the correlation strength r . As expected, we find that the asymptotic estimated 
covariance an is close to zero for r = 0 and it is close to the variance for r = 1. 
Figure 4(c) demonstrates that the interplay between covariances a\2 and variances 
an and #22 leads to a sublinear dependence of the asymptotic correlation coefficient 
p T = a 12/(^/^11^22) (Eq. (33)) in the large time bin limit, T = 25z s . 
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5 Conclusions 

Level-crossing phenomena occur in a variety of physical and biological sciences. 
In many of these situations coordination between level crossings of multiple cross- 
correlated Gaussian processes is of interest. Here, we focused on neuroscience and 
modeled the spikes of two cross-correlated neurons by two cross-correlated level- 
crossing processes. While crossings and extrema of one level-crossing process have 
been the focus of mathematical research, results describing the coordination of mul- 
tiple level-crossing processes are sparse and typically available only in specific and 
limited cases. Limits where level-crossing cross -correlations have been previously 
calculated are the weak and strong input correlation limit [3]. Here, we studied the 
case of two cross-correlated upward crossing processes and derived closed-form ex- 
pressions for their joint level-crossing coordination as well their joint count Gaus- 
sianity. Importantly, the results we present in this article are consistent with pre- 
viously reported limits but we now extended and generalized them. The two main 
results of our article are (1) closed-form explicit solution of the level-crossing cross- 
correlations and (2) the joint Gaussian limit of level-crossing counts. Our first result 
provides an explicit solution to v con d(^) = (si(t)s2(t + r))/ A /viV2 that is valid for 
all correlation strengths and which comprises previously obtained limits, see discus- 
sion in Sect. 3. The rate of level crossings by a one-dimensional Gaussian process 
is given by the prominent Rice's equation derived by Rice in the 1950s [11]. The 
solution we obtained for the level-crossing cross-correlation v C ond(^) extends the 
Rice rate to the joint rate of two correlated processes. Our second result proves the 
joint Gaussianity of level crossings for large bin sizes. The joint Gaussianity of spike 
counts is a highly desired property because if and only if two level-crossing counts 
are jointly Gaussian can zero count cross-correlation imply statistical independence. 
Notably, marginal Gaussianity of spike counts in each neuron combined with zero 
count cross-correlation is not sufficient to imply independence. Contrasting exam- 
ples of where X and Y variables are both marginally but not jointly Gaussian, have a 
zero cross-correlation but are not independent can be found in Sect. 5 in [8]. Count 
covariance and measures derived from it, such as the Pearson correlation coefficient, 
are computationally inexpensive and widely used as measures of statistical interde- 
pendencies [8]. Therefore, it is highly desirable to investigate the joint Gaussianity 
of level counts and thereby delimit the parameter space and mathematical conditions 
ensuring that independence can be implied from zero correlation coefficient. Notably, 
the joint Gaussianity of spike counts in bins of size T where T is much larger than 
the intrinsic time constant x s (T ^> r s ) also implies that models of multi-neuronal 
dynamics only need to consider the mean and variance of spike counts because all 
higher cumulants are zero. 

6 Materials and Methods 

6. 1 Proof of Proposition 3 . 1 

For simplicity of notation we adopt the following convention: (X\, Y\, X2, Y2) de- 
notes the vector (Vi(0), V/(0), V 2 (z), V^r)). In order to calculate (si(t)s 2 (t + r)), 
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defined in (14), we use the following regression model: 

Y\ =Ci\Xi +^2^2 + ^1, 

F 2 = ftXi+ftX 2 + 6 2 , 

where (£i , e 2 ) and (X\ , X 2 ) are independent. We use the covariance matrix in (4) and 
obtain 

27i4(T)27i3(T) -Z U (T)cr^ 



The conditional distribution C(Y\, Y 2 \Xi = \J/, X 2 = \J/) is a bivariate Gaussian dis- 
tribution 



f if{ai+a 2 ) ( cr €l (z) 2 Cov (€l , €2) (r) \ \ fi 
1 Wi+PiY \Cov (€u€2) (T) a €2 (r) 2 1 1 ' U V 

Using the regression system above, we write the conditional expectation 

E[ril{ FlG[ o,oo)}F2l{F 2 €[0,cx))}l^l =ir,X2 = f] 

= E [ 1 klG[-^(«l+«2),00]}l{62G[-^(^l+^2),00]} 

• (xl/ 2 (ai + a 2 )(P\ + ft) + ei^GBi + ft) + e 2 ^(ai + a 2 ) + *ie 2 )] 
= cr^r)^ (^)(^^E[1{Zig[«,oo)}1{Z 2 g[^,oo)}] - ^E[Zil{z ie [fl,oo)}l{Z2e[^oo)}] 
- aE[Z 2 l{z l e[a,oo)}^{Z 2 e[b,oo)}] + E[Zi l{ZiG[a,oo)}Z 2 l{z 2 e[£,oo)}]), (37) 
where Zi = -^U, Z 2 = a = -(^1+^2)) fe = -Wi+fe^ Now, we calcu- 

1 ffgj (r) ' z (7 e2 (r) ' a ei (r) ' o ei (r) 

late the four different expectations in Eq. (37) using Mehler's Formula (Lemma 10.7 
in [2]). First, we write 

E[Zil{ Zie[fl>00 )}Z2l{z 2 e[ft,oo)}] = y2 c n (a)c n (b)n ! ( C °^ (61 ' 62) ^ ) , (38) 

^ Vo- ei (T)or €2 (T)/ 

where c n (a) and c n (Z?) are the Hermite coefficients associated with the expectation 
E[Zi l{ZiG[a,oo)}^2l{Z 2 6[Z?,oo)}]- Then these Hermite coefficients are given by 

C( fl ) = ^1 f" z£L( e -*/*)6z. c n (b) = i=gl r z^(e~^)dz. 
n\V27t J a ^ n\V27t h &z n 

In particular, 

- for /? = 0, co(fl) = = 0(a), co(fe) = 
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- for n = 1, c\(a) = a</>(a) - 0(a) + 1, c\(b) = b(p(b) - 0(b) + 1, 

- for n > 2, c n (a) = ^( a H n _ x (a) + H n _ 2 (a)), c n (b) = ^(bH n _ x (b) + 
H n - 2 (b)). 

Analogously we obtain 

E[l {Zl €[ fl ,oo)}l{Z2€[^oo)}] = Yc n (a)c n (b)n\( , (39) 

^ \a ex (x)a ei (x) J 

with 

In particular, 

- for n = 0, c 0 (a) = l-cp(a) = 0(a), c 0 (b) = 0(b), 

- for n = 1, ci(a) = 0(a), ci(Z?) = 0(£>), 

- for n > 2, c„( fl ) = *$H n -i(a), c n (b) = *$H„-i(b). 

We also have 

E[l {Zl€ [ fl ,cx))}Z2l { z 2 €^oo)}] = y>(fl)cn(fr)*!( °^ (6l ' 62) ; rJ ) , (40) 

^ V<7 ei (r)or 62 (r)7 

with 

c„(a) = ^=r ^(e~ z2/2 )dz, c n (b) = r z^(e-^ 2 )dz. 

n\V2rc Ja dz n nl^/Tjr h &z n 

In particular, 

- lor n = 0, co(a) = 0(a), co(b) = 0(&), 

- for n = 1, ci(fl) = (j)(a), c x (b) = b</>(b) - 0(b) + 1, 

- for ft > 2, c n (fl) = i^Hn-iia), c n (b) = ^(bH n ^(b) + H n - 2 (b)). 

Finally, 

ElZil^g^oojl^g^ooj] = Vc w (a)c w (%!(^^) , (41) 

^ Vcr ei (rK 2 (r)/ 

with 

c„(a) = -<=gL rz^(e-^ 2 )dz, c n (b) = ^(e^dz. 
n\V2rc Ja dz n nl^/Tjr h &z n 

In particular, 

- lor n = 0, co(^) = (p(a), co(b) = 0(Z?), 

- for ft = 1, ci(tf) =a(p(a) - 0(a) + 1, c\(b) = (/)(b), 

- for ft > 2, c n (fl) = ^(aH n - X (a) + tf n _ 2 (fl)), Q (b) = ^H n ^(b). 
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Then, combining the four expectations (38)-(41), the associated Hermite coefficients, 
and Eq. (37), we obtain 

E[Fll{Fie[0,oo)}i2l{F 2 G[0,oo)}l^l =if,X2 = l/f] 

= (b0(b) -4>(b))(a0(a) -(P(a))a €l (x)a €2 {x) 
+ 0(a)0(b)CoY (6u€2) (t) 

+ E u C °] ( T 2) ^n i (4>(a)4>(b)H n -2(a)H n - 2 (b)). (42) 

n>2 1 z 

Note that the first two terms in Eq. (42) correspond to orders n = 0 and n = 1 , re- 
spectively. Denoting E[Fi 1{Fig[0,oo)}^21{F 2 g[0,oo)}I^1 = f,X 2 = f] = C {a , b) (z) we 
find {s\(t)s 2 (t + ?)) = C( a ,b)(^)Pz(^, VO- Here, C{ a ,b)(^) is a uniformly convergent 
series. □ 

6.2 Zero Time Lag Correlations 

In this section we derive the zero lag spike correlation using the Gaussian probabil- 
ity integrals in Eq. (15). Following the previously introduced notation the spiking 
threshold level in a neurons is ^ and variance ay i we can write 

VcondM = {Sl(t)s 2 (t))/(VIV 2 ) 
j poo poo 

= -== / / Vi • m + T) P rW, V, Vi, y 2 )d\>idy 2 (r). 

Now, we substitute the variables: 



E = 



r = 



y j2ay l ay 2 (ray l ay 2 -\-ay 1 ay 2 ) 
joy 2 lQy x + y/cr^/cr^Vjt + r) 

V2 Jo^y 2 ~ rOy x Oy 2 c"(x) 

A _ VKay 2 -QTVi) 

yl2ay l Gy 1 (py l Gy 1 — rOy x Oy 2 ) 

\f2loy x Gy 2 + ror yi ory 2 c" (r ) 

The correlation matrix ^ ^ ^ for r = 0 is a four-dimensional identity matrix. 

2 2 

Vcond(O) = 



J.ooL F V 2(1 +r) 



)i^i7(i-'- 2 ) 
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+ ^[(l+r)-A 2 (l-r)]Erfc(^ 

if 2 ((T V i +ov 2 ) 2 



VVOTT) 



• exp 



4a y l Oy 2 {rOy l Oy 2 + Oy x Oy 2 ) 
ir 2 (a V2 -a Vi ) 2 A 2 



2 " vu - 'dA (43) 



4<Ty l <Ty 2 ((7y l <7y 2 — YOy x Gy 2 ) 2 

Solving this integral we obtain Eq. (23). 
6.3 Proof of Theorem 4.1 

To prove joint Gaussianity, it is instructive to briefly recapitulate the one-dimensional 
Central Limit Theorem. We start by writing 

oo oo r j 

U^ T] (xlfi) = J2J2 d T^ i)ak / ^(Xtis^H^X^ds, (44) 

j=0k=0 ^° 

with dj^iijfi) = jj(j)(\Jfi)Hj(\lfi) and a k = ^ / 0 °° xHk(x)(j)(x)dx. Defining the level 
count deviation for neuron i by S l (T) we obtain 

= jt~7= [ T df(f i )a k H j (X i (s))H k (X f i (s))ds (45) 
q=l ^TJo k+j=q 

= £4f f T G^XtisXX'^ds (46) 



= ^^(r.Xi.x;), (47) 

where G l q (x\,X2) = ^j i+ j =q dJ\\lri)a] i Hj(x\)H] i (x2). A Gaussian distribution is a 
stable limit distribution for a sum of independent finite variance variables. Therefore, 
all that is left to prove is that contributions q ^ q' are independent and have finite 
variance. From Mehler's Formula we recognize that the contributions for q ^ q' are 
independent. The finite variance follows from the observation that for all q the vari- 
ance of G l q (Xi (s), X\ (^)) is proportional to the expectation of a product of four Her- 
mite polynomials, which has been proven to be finite (Theorem 10.10 in [2]) if the 
conditions of Theorem 4.1 are satisfied. 

Now, we address the joint Gaussianity. A random vector is jointly Gaussian if and 
only if any linear combination of its components has a univariate normal distribution 
(see [4, 32, 33]). Thus, we need to prove that the sequence a\S l (T) + ct2S 2 (T) is 
asymptotically Gaussian and satisfies a Central Limit Theorem, for all a\ , ai e R. We 
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start from Eq. (45) and consider a truncated series for S 1 q (T) = J2q=\ Jq (^» ^1 » ^i) 
consisting of the first Q terms and denote by || • || the norm in L 2 {£2). First, we show 
that the remainder is bounded 



lim 

T—^oo 



i=l 



i=l 



< m\ 



crl x {q) + \a 2 \ 



E °-l 2 («)- 

M <z=0+i 



(48) 



As 2 grows, this difference diminishes such that if a\S l Q(T) + ci 2 S 2 q(T) is Gaussian 
for large Q then this will imply that a\S l (T) + a 2 S 2 (T) is Gaussian, too. Now, we 
only need to show Gaussianity of a\S l Q(T) + ci 2 S 2 q(T). Using a modified version of 
the Breuer-Major Theorem (Sect. 6.4), we know that for each q 



t l jl(T,X l ,X> l )+a 2 J 2 q {T,X 2 ,X' 2 ) iV(0,a? ) 



where a? is given in Eq. (51). The same theorem implies that a\J^(T, X\, X[) + 

a 2 J 2 (T, X 2 , X f 2 ) and a\J^ r {T, X\, X[) + a 2 J 2 f (T, X 2 , X f 2 ) are asymptotically inde- 
pendent if q ^ q' . Thus we obtain for any Q > 1 

^ Q (T)+a 2 Sl(T)^Af(o,j24 q )- 
To calculate the count covariances in Eq. (26) we use Lemma 10.7 in [2] and obtain 

oo q q 

ai J = 2 EE E 4 l \(iri)d^(irj)a k a^ 

q=lk=0 k'=0 
poo 

x / E[^_^(ZK0))^(z;(0))^_^(Zy(5))//^(x;.(5))]d^ (49) 

Applying the Mehler Formula to a four-dimensional Gaussian random vector 
(Lemma 10.7 in [2]) we find 

E[H q - k (Xi (0)) H k (X[ (0)) H q ,_ v (Xj (s)) H t (X) (s))] 
[0, 



if q y^q\ 

di\d 2 \d 3 W 

•c(s) dl c f (s) d2 (-c f (s)) d i(-c"(s)) d \ if q = q f . 



V ) J 2^(di,d 2 ,d 3 ,d 4 )eZ 



Here, 8ij is the Kronecker delta, Z is the set of di 's satisfying: di >0,d\-\-d 2 = q—k, 
di + d4 = k, d\ + di = q — k\ and d 2 + J4 = k' . We thus can express the covariances 
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of the spike counts as 



q=l k=0k'=0 



poo 



(q - k)lk\(q - k')\k' 



d x \d 2 \d 3 \d A \ 
.c(s) d ic\s) d2 (-c\s)) d3 (-c"(s)) d4 ds. 
This is the result reported in Theorem 4.1. □ 

6.4 Modified Breuer-Major Theorem 

Here, we adapt the Breuer-Major Theorem [34] to show that the bivariate vector 
(Jq(T, X u X[), J%(T, X 2 , X f 2 )) is Gaussian. 

Theorem 6.1 We consider two zero mean and unit variance Gaussian processes 
Xf(t), which are described by the properties in Sect. 2 and a correlation func- 
tion E[Xi(0)Xj(t)] = Cij(t), where i, j e {1, 2}. For all functions Gj(-, •) that fulfill 
E[G/(X/(0), X ■(()))] = 0 and E[G?(X f (0), X ■(()))] < oo and two real constants a { 
the following integral is Gaussian and we have 



where 



l — [ T (a 1 G l (X 1 (t),X[(t)) + a 2 G 2 (X 2 (t),X f 2 (t)))dt -^U tf(0,cr£), 

'T JO T^oo 



poo 

al = 2aj / E[Gi(Xi(0),Zi(0))Gi(Zi(f),x;(f))]df 
Jo 

poo 

+ 2a\ \ E[G 2 (X 2 (0), X' 2 (0))G 2 (X 2 (t), X' 2 (t))]dt 
Jo 

poo 

■4a x a 2 / E[Gi(X 1 (0),X / 1 (0))G 2 (X 2 (0,4(0)]dL (51) 
Jo 



Proo/ We start by considering the Hermite expansion of the functions G;(-, •) and 
write 

G 

Gi(Xi(f),^(f))= lim Y Y Ck^Hk^X^Hk^X^t)) (52) 

Q^-oo*—^ 

q=l ki+k2=q 

= lim Gp, (53) 
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where Gp denotes the sum over q in Eq. (52) truncated at Q, and the convergence is 
in L 2 (£2). We write 

[ T a 1 G 1 (X 1 (t),X[(t))+a 2 G 2 (X 2 (t),X f 2 (t))dt (54) 
= lim -L [ T (a 1 G¥(X 1 (!) 9 X' 1 (t))+a 2 G$ (X 2 (t), X f 2 (t)))dt. (55) 

<2^oo Jo 

Now, to prove the Gaussianity of Eq. (54) it is sufficient to prove the asymptotic 
Gaussianity of (aiGp (Xi(t), X[(t)) + a 2 G$(X 2 (t), X' 2 (t))). This is sufficient be- 
cause 



Kq(T) = 



j r 2 

v ^ Jo 

^Y>\-jfh [(Gi(Xi(t),X' { (t)) -G?(Xi(t),X' { (t)))] 
where || • || is the norm in L 2 {£2). Each of the terms is bounded, 

II T 2 

f [Gi(Xi(t), Xfr)) - G?(Xi(t), X' t (t))]dt 

II vr Jo 

^Ef(i-f) E Wife" 
+ X)f O-^W E ^ 2 , Gi w 

G + l " V 7 £l+£ 2 =4 

oo 

<«E E c w 2 ,g,w 

G+l&l+&2=? 

+ Ef( 1 -^V )G E 4 lk2 , Gi w> 



(56) 



(57) 



(58) 



(59) 



(60) 



where ^(f) = max(|cn(OI + |ci2(0l> k2i(0l + 1^22(01) an d ^ is a positive real 
number such that < 1 whenever t > a. Since Eq. (60) is a vanishing series, 
limg^oo lim^^oo Kq(T) = 0. Now, we address the Gaussianity of (a\Gf (X\(t), 
X[(t))+a 2 G%(X 2 (t), X r 2 (t))). A result of Kratz and Leon, Lemma 3, p. 653 in [35], 
implies that 



L i (T) = -^j*G?(X i (t),X' i (t))dt 
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can be approached in L 2 (Q) by the sequence of € -dependent processes Xf(t) 
L^ S (T) = -j= / 0 r G? (Xf(0, (Xf) 7 (0)dr. The ^-dependent processes for i = 1,2, 
are defined as 

/oo 
^(/v * AW) 7 (VT^7dW/W + V^dW 3 (A)), (61) 
-oo 

where * denotes a convolution, /J e is /3 e (t) = j/Kf)> being a positive function 
with |/3(A)dA < oo, j = 1, 2 and such that its Fourier transform has support 

in [—1, 1], where E[X?(0)X?(r)] = Cjj(r)P(sr) these processes are ^-dependent. 
Then 

II T 2 

II v r Jo i=1 

2 

<J2\^i\\W(T)-L^(T)l (62) 



2 



Urn Urn y^\oii\ II V (T) - L i e (T) II = 0. (63) 



i = l 



This follows from the Central Limit Theorem for ^ -dependent random vectors and 
concludes the proof. □ 
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